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Abstract 

The Generalized Fluid System Simulation Program (GFSSP) is a finite-volume based 
general-purpose computer program for analyzing steady state and time-dependent flow 
rates, pressures, temperatures, and concentrations in a complex flow network. The program 
is capable of modeling real fluids with phase changes, compressibility, mixture 
thermodynamics, conjugate heat transfer between solid and fluid, fluid transients, pumps, 
compressors, flow control valves and external body forces such as gravity and centrifugal. 
The thermo-fluid system to be analyzed is discretized into nodes, branches, and conductors. 
The scalar properties such as pressure, temperature, and concentrations are calculated at 
nodes. Mass flow rates and heat transfer rates are computed in branches and conductors. 
The graphical user interface allows users to build their models using the ‘point, drag, and 
click’ method; the users can also run their models and post-process the results in the same 
environment. The integrated fluid library supplies thermodynamic and thermo-physical 
properties of 36 fluids, and 24 different resistance/source options are provided for modeling 
momentum sources or sinks in the branches. Users can introduce new physics, non-linear 
and time-dependent boundary conditions through user-subroutine. 

1. Introduction 

The need for a generalized computer program for thermo-fluid analysis in a flow network 
has been felt for a long time in Aerospace Industries. Designers of thermo-fluid systems 
often need to know pressures, temperatures, flowrates, concentrations, and heat transfer 
rates at different parts of a flow circuit for steady state or transient conditions. Such 
applications occur in propulsion systems for tank pressurization, internal flow analysis of 
rocket engine turbo-pumps, chilldown of cryogenic tanks and transfer lines and many 
other applications of gas-liquid systems involving fluid transients and conjugate heat and 
mass transfer. Computer resource requirements to perform time-dependent three- 
dimensional Navier-Stokes Computational Fluid Dynamic (CFD) analysis of such 
systems are prohibitive and therefore are not practical. A possible recourse is to 
construct a fluid network consisting of a group of flow branches such as pipes and ducts 
that are joined together at a number of nodes. They can range from simple systems 
consisting of a few nodes and branches to very complex networks containing many flow 
branches simulating valves, orifices, bends, pumps and turbines. In the analysis of 
existing or proposed networks, node pressures, temperatures and concentrations at the 
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system boundaries are usually known. The problem is to determine all internal nodal 
pressures, temperatures, concentrations and branch flow rates. Such schemes are known 
as Network Flow Analysis methods and they use largely empirical information to model 
fluid friction and heat transfer. 

The oldest method for systematically solving a problem consisting of steady flow in a 
pipe network is the Hardy Cross method [1], The Hardy Cross method works well for 
hand calculation but experiences slow convergence for large circuit, The network 
analysis method has been widely used in thermal analysis codes (SINDA/G [2] and 
SINDA/FLUINT [3]) using an electric analog. The partial differential equation of heat 
conduction is discretized into finite difference form expressing temperature of a node in 
terms of temperatures of neighboring nodes and ambient nodes. The set of finite 
difference equations are solved to calculate temperature of the solid nodes and heat fluxes 
between the nodes. GFSSP [4] uses a “pressure -based” finite volume method [5] as the 
foundation of its numerical scheme. 

This paper provides a brief overview of GFSSP’s network definition, data structure, 
mathematical formulation, thermodynamic property program and program structure and 
describes additional capabilities of version 6. The paper also describes several validation 
and verification effort. 

2. GFSSP Overview 

2.1 Network Definition 


A fluid network using fluid and solid nodes. The fluid circuit is constructed with boundary 
nodes, internal nodes and branches (Figure 1) while the solid circuit is constructed with 
solid nodes, ambient nodes and conductors. The solid and fluid nodes are connected with 
solid-fluid conductors. Users must specify conditions, such as pressure, temperature and 
concentration of species at the boundary nodes. These variables are calculated at the 
internal nodes by solving conservation equations of mass, energy and species in 
conjunction with the thermodynamic equation of state. Each internal node is a control 
volume where there are inflow and outflow of mass, energy and species at the boundaries 
of the control volume. The internal node also has resident mass, energy and concentration. 
The momentum conservation equation is expressed in flowrates and is solved in branches. 
At the solid node, the energy conservation equation for solid is solved to compute 
temperature of the solid node. Figure 1 shows a schematic and GFSSP flow circuit of a 
counter- flow heat exchanger. Hot nitrogen gas is flowing through a pipe, colder nitrogen 
is flowing counter to the hot stream in the annulus pipe and heat transfer occurs through 
metal tubes. The problem considered is to calculate flowrates and temperature distributions 
in both streams. 
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Figure 1 - A Typical Flow Network consists of Fluid Node, Solid Node, Flow 

Branches and Conductors 


2.2 Data Structure 

GFSSP has a unique data structure (Figure 2) that allows constructing all possible 
arrangements of a flow network with no limit on the number of elements. The elements of 
a flow network are boundary nodes, internal nodes and branches. For conjugate heat 
transfer problems, there are three additional elements, solid node, ambient node and 
conductor. The relationship between a fluid node and a branch as well as a solid node and 
conductor is defined by a set of relational geometric properties. For example the relational 
geometric properties of a node are number and name of branches connected to it. With the 
help of these properties, it is possible to define any structure of the network as it progresses 
through every junction of the network. The positive or negative flow direction is also 
defined locally. Unlike structured co-ordinate system, there is no global definition of flow 
direction and origin. The development of a flow network can start from any point and can 
proceed in any direction. 
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Figure 2 - Data Structure of the Fluid-Solid Network has Six Major Elements 

2.2 Mathematical Formulation 

GFSSP solves the conservation equations of mass and momentum in internal nodes and 
branches to calculate fluid properties. It also solves for energy conservation equations to 
calculate temperatures of solid nodes. Table 1 shows the mathematical closure that 
describes the unknown variables and the available equations to solve the variables. 
Pressure, temperature, species concentration and resident mass in a control volume are 
calculated at the internal nodes whereas the flowrate is calculated at the branch. The 
equations are coupled and non-linear. Therefore, they are solved by an iterative numerical 
scheme. GFSSP employs a unique numerical scheme known as Simultaneous Adjustment 
with Successive Substitution (SASS) which is a combination of Newton-Raphson and 
Successive Substitution methods. The mass and momentum conservation equations and 
the equation of state are solved by the Newton-Raphson method while the conservation of 
energy and species are solved by the successive substitution method. 
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Table 1 Mathematical Closure 


Unknown Variables 

1. Pressure 

2. Flowrate 

3. Fluid Temperature 

4. Solid Temperature 

5. Species Concentrations 

6. Fluid Mass (Unsteady Flow) 
2.3 Fluid Properties 


Available Equations to Solve 

1. Mass Conservation Equation 

2. Momentum Conservation Equation 

3. Energy Conservation Equation of Fluid 

4. Energy Conservation Equation of Solid 

5. Conservation Equations for Species 

6. Thermodynamic Equation of State 


GFSSP is linked with two thermodynamic property programs, GASP & WASP [6, 7] and 
GASPAK [8] that provide thermodynamic and thermo-physical properties of selected 
fluids. Both programs cover a range of pressure and temperature that allows fluid 
properties to be evaluated for liquid, liquid-vapor (saturation) and vapor region. GASP 
and WASP provide properties of twelve fluids (Table 2). GASPAK includes a library of 
thirty six fluids (Table 3). 


Table 2 - Fluids Available in GASP & WASP 


Index 

Fluid 

Index 

Fluid 

1 

HELIUM 

7 

ARGON 

2 

METHANE 

8 

CARBON DIOXIDE 

3 

NEON 

9 

FLUORINE 

4 

NITROGEN 

10 

HYDROGEN 

5 

CARBON MONOXIDE 

11 

WATER 

6 

OXYGEN 

12 

RP-1 
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Table 3 - Fluids Available in GASPAK 


Index 

Fluid 

Index 

Fluid 

1 

HELIUM 

19 

KRYPTON 

2 

METHANE 

20 

PROPANE 

3 

NEON 

21 

XENON 

4 

NITROGEN 

22 

R-11 

5 

CO 

23 

R12 

6 

OXYGEN 

24 

R22 

7 

ARGON 

25 

R32 

8 

CG 

26 

R123 

9 

PARAHYDROGEN 

27 

R124 

10 

HYDROGEN 

28 

R125 

11 

WATER 

29 

R134A 

12 

RP-1 

30 

R152A 

13 

ISOBJTANE 

31 

NITROGEN TRIFLUORIDE 

14 

BUTANE 

32 

AMMONIA 

15 

DEUTERIUM 

33 

IDEAL GAS 

16 

ETHANE 

34 

HYDROGEN PEROXIDE 

17 

ETHYLENE 

35 

AIR 

18 

HYDROGEN SULFIDE 




2.4 Program Structure 

GFSSP has three major parts (Figure 3). The first part is the Graphical User Interface, 
VTASC (Visual Thermofluid Analyzer of Systems and Components). VTASC allows 
users to create a flow circuit by a ‘point and click’ paradigm. It creates the GFSSP input 
file after the completion of the model building process. It can also create a customized 
GFSSP executable by compiling and linking User Subroutines with the solver module of 
the code. Users can run GFSSP from VTASC and post process the results in the same 
environment. The second major part of the program is the Solver and Property Module. 
This is the heart of the program that reads the input data file, generates the required 
conservation equations for all internal nodes and branches with the help of thermodynamic 
property data. It also interfaces with User Subroutines to receive any specific inputs from 
users. Finally, it creates output files for VTASC to read and display results. The User 
Subroutine is the third major part of the program. This consists of several blank subroutines 
that are called by the Solver Module. These subroutines allow the users to incorporate any 
new physical model, resistance option, fluid etc. in the model. 
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Figure 3 - GFSSP’s Program Structure showing the interaction of three major 

modules 


2.5 Resistance Option 

In network flow analysis code, flow resistances are modeled by empirical laws. These 
empirical laws have been incorporated to model flow resistances for pipe flow, orifices, 
valves and various pipe fittings. GFSSP models these flow resistances in the momentum 
conservation equation as a friction term. There are twenty- four different resistance options 
available to users to choose from. There is also a provision for introducing a new resistance 
option through user subroutines. The available resistance options are shown in Table 4. 

2.6 Graphical User Interface 

GFSSP’s Graphical User Interface provides the users a platform to build and run their 
models. It also allows post-processing of results. The network flow circuit is first built 
using three basic elements: boundary node, internal node and branch. Then the properties 
of the individual elements are assigned. Users are also required to define global options of 
the model that includes input/output files, fluid specification and any special options such 
as rotation, heat exchanger etc. During execution of the program, a run manager window 
opens up and users can monitor the progress of the numerical solution. On the completion 
of the run, it allows users to visualize the results in tabular form for steady-state solutions 
and in graphical form for unsteady solutions. It also provides an interface to activate and 
import data to the plotting program, WINPLOT [9] for post processing. 
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Table 4 - Resistance Options in GFSSP 


Optio 

n 

Type of 
Resistance 

Input 

Parameters 

Opti 

on 

Type of 
Resistance 

Input Parameters 

1 

Pipe Flow 

L (in), D (in), 

s/D 

13 

Common Fittings 
and Valves (Two 
K Method) 

D (in), Ki, K 2 

2 

Flow Through 
Restriction 

Cl, A (in 2 ) 

14 

Pump 

Characteristics 1 

Ao, Bo, Co, A (in 2 ) 

3 

Non-circular Duct 

a (in), b (in) 

15 

Pump Power 

P (hp), p, A (in 2 ) 

4 

Pipe with 
Entrance and Exit 
Loss 

L (in), D (in), 

8/D, Ki, Ke 

16 

Valve with Given 
Cv 

Cv , A 

5 

Thin, Sharp 
Orifice 

Di (in), D 2 (in) 

17 

Joule-Thompson 

Device 

Lohm, Vf, kv, A 

6 

Thick orifice 

L (in), Di (in), 
D 2 (in) 

18 

Control Valve 

See Example 12 data 
file 

7 

Square Reduction 

Di (in), D 2 (in) 

19 

User Defined 

A (in 2 ) 

8 

Square Expansion 

Di (in), D 2 (in) 

20 

Heat Exchanger 
Core 

Af (in 2 ), A s (in 2 ), A c 
(in 2 ), L(in), K c , K e 

9 

Rotating Annular 
Duct 

L (in), r 0 (in), 
r; (in), N (rpm) 

21 

Parallel Tube 

L (in), D (in), s/D, n 

10 

Rotating Radial 
Duct 

L (in), D (in), 
N (rpm) 

22 

Compressible 

Orifice 

Cl, A (in 2 ) 

11 

Labyrinth Seal 

n (in), c (in), m 
(in), n, a 

23 

Labyrinth Seal, 
Egli Correlation 

r; (in), c (in), m (in), n, 

a 

12 

Flow Between 
Parallel Plates 

n (in), c (in), 
L (in) 

24 

Fixed Flow 

Flow (lbm/s), A (m 2 ) 


2 


1 Pump characteristics are expressed as Ap = Ao + Bo m + Co m 

A p - Pressure rise, lbf/ft 2 

m - Flow rate, lbm/sec 


2.7 Example Problems 

Several example problems have been included to aid users to become familiar with 
different options of the code. The example problems also provide the verification and 
validation of the code by comparing code’s predictions with analytical solution and 
experimental data. The example problems are: 
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(1) Simulation of a flow system consisting of a pump, valve, and pipe line. 

(2) Simulation of a water distribution network. 

(3) Simulation of compressible flow' in a converging-diverging nozzle. 

(4) Simulation of the mixing of combustion gases and a cold gas stream. 

(5) Simulation of a flow system involving a heat exchanger. 

(6) Radial flow on a rotating radial disk. 

(7) Flow' in a long-bearing squeeze film damper. 

(8) Simulation of the blowdown of a pressurized tank. 

(9) A reciprocating piston-cylinder. 

(10) Pressurization of a propellant tank. 

(11) Power balancing of a turbopump assembly. 

(12) Helium pressurization of LOX and RP-1 propellant tanks. 

(1 3) Steady state and transient conduction through a circular rod, with convection. 

(14) Chilldown of a short cryogenic pipe line. 

(1 5) Simulation of fluid transient following sudden valve closure. 

(16) Simulation of pressure regulator downstream of a pressurized tank. 

(1 7) Simulation of flow regulator downstream of a pressurized tank. 

(18) Subsonic Fanno flow. 

(19) Subsonic Rayleigh flow. 

(20) Modeling of closed cycle liquid metal (lithium) loop with heat exchanger to heat helium 
gas, 

(21 ) Internal flow in a turbopump. 

(22) Simulation of a fluid network with fixed flow rate option. 

(23) Helium-assisted, buoyancy -driven flow in a vertical pipe carrying LOX with ambient 
heat leak. 

(24) Simulation of relief valve in a pressurized tank. 

(25) Two-dimensional recirculating flow in a driven cavity. 

(26) Fluid transients in pipes due to sudden opening of valve. 

(27) Boiling water reactor. 

(28) No-vent tank chill and fill model. 

(29) Self-pressurization of a cryogenic propellant tank due to boil-off. 

(30) Modeling solid propellant ballistic with GFSSP, 

The Tables (5a & 5b) show the particular features of each example problem. For 
example, “Conjugate Heat Transfer” option has been used in Examples 13, 14, 23, 26 and 
29. 
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Table 5a 


Table 14a. Use of various options in example problems — examples 1-15. 


Feature 

Example 

1 

2 

3 

4 

S 

G 

7 

3 

9 

10 

11 

12 

13 

14 

IS 

Conjugate heat transfer 













13 

14 


Constant properly 


2 





7 









Cydic boundary 
















Rxed mass flow 
















Row regulator 
















Gravity 

1 















Heat exchanger 





5 






11 





Ideal gas 








a 








Long inertia 



3 



6 






12 




Ruid mixture 




4 






ID 


12 




Model import 
















Moving boundary 







7 


9 







Multilayer insulation 
















Multidimensional flow 
















Nondrajlar duel 







7 









Phase change 














14 


Pressurization (tank) 










ID 


12 




Pressure regulator 
















Pressure relief valve 
















Pump 

1 











12 




Solid rocket motor 
















Tuftopump 











11 





Turbopfmp-Hntemal flew 
















Unsteady 








8 

9 

ID 


12 


14 

15 

User fluid 
















User subroutine 










ID 


12 




VafveO 10 















15 

Variable geometry 









9 







Ruid transient 
(water hammer) 















15 
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Table 5b 


Table 14b. Use of various options in example problems — examples 16-30. 


Feature 

Example 

16 

17 

IS 

1$ 

20 

21 

22 

23 

24 

25 

26 

27 

28 

29 

30 

Conjugate heat transfer 








23 





28 

29 


Constant property 










25 






Cydic boundary 





20 











Fixed mass flow 







22 






28 



Row regulator 


17 














Gravity 








23 




27 

28 

29 


Heat exchanger 





20 











Ideal gas 

18 

17 













30 

Long media 



18 

19 








27 



30 

Ruid mixtire 








23 








Model import 








23 








Moving boundary 











28* 





Multilayer insulation 














29 


Multidimensional flow 










25 






Norrarcular dud 
















Phase change 












27 

28 

29 


Pressurization (tank) 














29 


Pressure regulator 

18 















Pressure relief vafve 









24 







Pump 
















Solid rocket motor 















30 

Turhopunp 
















Turbopimp — internal flow 






21 










Unsteady 

16 

17 





22 

23 

24 


20 


28 

29 

30 

User fluid 





20 







27 




User subroutine 



18 

19 

20 






26 


28 

29 

30 

VaJveO/C 











26 





Variable geometry 











26* 





Ruid transient 
(water hammer) 











26 






•Variable geometry and moving boundary handled by Use* Subrowme 
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3. 


Additional Capabilities of Version 6 


The additional capabilities of Version 6 include improved modeling of fluid mixtures, 
multi-dimensional flow capability in a system level flow network model, extension of the 
thermo-dynamic property package, extension of the pressure and flow regulator option, 
and inclusion of the relief valve and fixed flowrate options. 

3.1 Mixture modeling 

The mixture modeling capability in earlier version of GFSSP did not allow phase change 
for any component of the mixture. This limitation was due to the fact that the energy 
conservation equation was expressed as a product of specific heat and temperature instead 
of enthalpy. In many liquid propulsion applications, phase change in mixture is a common 
occurrence when cryogenic propellants mix with inert gas such as helium or nitrogen. This 
limitation was removed in Version 6 by introducing an additional option where the energy 
equation of each species was expressed in terms of enthalpy. This option allows change 
of phase for any component of the mixture. 

3.2 Multi-dimensional flow capability 

Version 6 allows the user to model multi-dimensional flow in a fluid network. Network 
modeling usually applies to a system where a one-dimensional momentum equation is 
sufficient to characterize the flow. However, in some situations, such as a stratified 
cryogenic tank, the one-dimensional flow assumption is not realistic. Multi-dimensional 
flow modeling in a system level code will eliminate the need to integrate with a CFD code 
which has not yet been proved to be a practical solution. Multi-dimensional flow modeling 
in a system level code is a viable alternative to address such a need. Multi-dimensional 
capability in GFSSP has been verified by comparing its predictions with classical 
numerical fluid dynamics problems such as flow in a driven cavity. 


3.3 Improvements and Extension in Thermodynamic Property Routines 

There has been significant improvement and extension of the thermodynamic property 
routines in Version 6. The thermodynamic property routines have been rewritten to 
introduce universal property call subroutines. In this process the fluid property code has 
been reduced by nearly 1000 lines. With the introduction of the new property call routines, 
user subroutines can make property calls during any stage of the computation. 

The user supplied fluid property table has been extended to include an optional saturation 
table. It may be noted that earlier version of GFSSP did not have the capability to model 
phase change with user supplied fluid property tables. Utility programs have been 
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developed to generate user supplied thermodynamic property tables from tables generated 
by REFPROP [10] and are included in the GFSSP installation package. 


3.4 Extension of pressure and flow regulator option 

A marching algorithm [11] for modeling pressure and flow regulator has been introduced 
in Version 6. The earlier option of iterative algorithm is still available. The marching 
algorithm is more economic and allows the use of multiple regulators in a given flow 
circuit. 


3.5 Inclusion of relief valve and fixed flowrate option 

Version 6 has also the capability of modeling relief valve and a fixed flowrate option, ft 
may be recalled that GFSSP’s mathematical formulation requires pressures to be specified 
at boundary nodes and flowrates are calculated by solving the momentum conservation 
equation. In this version the user can specify a given flowrate in a branch connected with 
a boundary node. This option is available for both steady and unsteady flow. 


4. Validation & Verification 

Several code validation efforts were completed during the development of Version 6. The 
applications considered for code validation include fluid transient, chilldown of cryogenic 
transfer line, self-pressurization of cryogenic tank, and no vent chill and fill of cryogenic 
tank. 

4. 1 Two-dimensional recirculating flow in a square cavity 

In this example, two-dimensional recirculating flow in a square cavity [12] has been modeled 
using GFSSP’s multi-dimensional flow calculation capability. In a square cavity the flow is 
induced by shear interaction at the top wall as shown in Figure 4. The length of each wall 
is 12 inches. The density of the fluid is assumed constant at 1.00 lbm/ft 3 , and the viscosity 
of the fluid is assumed to be 1.00 lb m /(ft sec). The bottom and side walls are fixed. The 
top wall is moving to the right at constant speed of 100 ft/sec. The corresponding Reynolds 
number for this situation is Re = 100. 

The GFSSP model (Figure 5) of the driven cavity consists of 50 nodes (49 of which are 
internal) and 84 branches. System model is shown in Figure 5a. The expanded component 
(square cavity) is shown in Figure 5b. 
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Figure 4- Flow in a Shear Driven Square Cavity 
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(a) System Network with expandable grid (Element 1) 


14 




(b) Expanded two-dimensional Cartesian grid of Element 1 


Figure 5- Two-dimensional Cartesian grid generation in VTASC 

Figure 6 shows a comparison between the benchmark numerical solution and GFSSP 7x7 
node model velocity profiles along a vertical plane at the horizontal midpoint. As can be 
seen in Figure 6, the results of this crude GFSSP model compare very favorably with the 
benchmark numerical solution of Burggraf [13], 
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Figure 6. Shear Driven Square Cavity Centerline Velocity Distribution 

The predicted velocity field and pressure contours are shown in Figure 7. The 
recirculating flow pattern and stagnation of flow near the top right comer are clearly 
shown in the figure. The predicted stream traces from calculated velocity field is shown 
in Figure 8. 


Frame 001 1 16 Dec 2010 | 



Figure 7 - Predicted Velocity Field and Pressure Contours 
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Figure 8 - Predicted Stream Traces in the Driven Cavity 


4.2 Fluid Transient 

This example deals with water hammer in a pipe with entrapped air. GFSSP results are 
validated against experimental data available in the literature [14]. A long pipe of 
diameter 1 .025 inches is attached to a reservoir of water at one end and closed at the other 
end with some entrapped air in the other end. A ball valve separates the water from the 
air as shown in Figure 9 below. The ball valve is closed until about 0.15 second, and then 
gradually opens to 100% at about 0.4 second. This example has been set up according to 
the experimental study done by Lee et al[14]. The two most important controlling 
parameters for this problem are the reservoir pressure (pr) and the fractional air length 
present in the pipe as compared to the total pipe length (a g = Lg/Li ). The initial length 
for the water volume in the pipe (Li) is fixed to 20 ft, and initial length of air column in 
the pipe (L g ) varies from a low of 1.23 ft to 16.23 ft, the value of a ranging from 0.0579 
to 0.448 respectively. The ratio of reservoir pressure to the initial pressure of the 
entrapped air (Pr = pR/p a tm) varies in the range of 2 to 7, i.e. the reservoir pressure (p R ) 

range being 29.4 psi to 102.9 psi. The objective of this study is to predict the transient 
pressure at different points along the length of the pipe. 
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Figure 229. Schematic of the Water Pipe with Entrapped Air 


The GFSSP model to represent the flow of water in the pipe is shown in Figure 10. The 
pipe sector of length 20 feet (only the water column) is divided into ten uniform pipe 
segments and one restriction separating twelve nodes. Boundary Node 1 represents the 
tank (reservoir). A user subroutine interfaces Node 12 to an unseen pseudo control 
volume containing air only. The pseudo control volume has a fixed mass of air, but the 
volume changes as it is pressurized owing to the fluctuation of pressure at node 12. 
Thereby the volume of node 12 changes as the volume of the imaginary control volume 
changes. The volume change in node 12 is computed by a volume balance between the 
volume of water and the volume of the entrapped air. The total volume (Vtot = Vair+ V12) 
remains constant, and must be equal to the initial total volume (since the pipe is closed at 
the other end). 


H— «-HZI— ■HZI — ■-HI} 

12 23 34 45 

VQC 



■— H—— S 

1112 1011 910 

Figure 10. GFSSP Model of Sudden Valve Opening Experiment of Lee & Martin 54 


For the numerical solution a time step of 0.01 s has been used. The operating conditions 
are: PR = 7, and a g = 0.45). Figure 1 1 compares GFSSP’s predicted pressure at node 12 
with that of the experimental data points of Lee and Martin. The predicted results 
compared very well, and even though the peak pressure amplitude differs by about 7%, 
the frequencies of pressure oscillations matched very well. 
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Figure 11. Comparison of GFSSP and Experimental Data 

A Fast Fourier Transform (Figure 12) has been conducted in the numerical model to 
predict the different modal frequencies of the pressure transient and also compared with 
the experimental data. More details of this example problem is available in Reference 
[ 15 ]. 



Figure 12. Fast Fourier Transform for Modal Frequencies 
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4.3 Chilldown of Cryogenic Transfer Line 

For this example, the chilldown of cryogenic pipeline to validate GFSSP’s transient 
conjugate heat transfer capability has been selected. In the 1960s the National Bureau of 
Standards (NBS) conducted a series of chilldown experiments on a cryogenic transfer 
line [16]. The test setup (Figure 13) is a vacuum-jacketed 200 foot long copper pipe of 
5/8 inch inner diameter. A pressurized 80 gallon dewar feeds liquid hydrogen into the 
pipe that is initially at ambient temperature. The wall temperature is measured at four 
thermocouple stations at distances of 20, 80, 141, and 198 feet from the inlet. 

When the fluid touches the relatively warm pipe walls, heat transfer causes the liquid 
cryogen to boil and the pipe wall temperature to decrease. Eventually the pipe chills 
down to the liquid temperature, and the liquid front gradually travels further down the 
pipeline. At the outlet of the pipeline, vapor exits to the atmosphere. 

The NBS experiments were conducted with liquid hydrogen (LH2) and liquid nitrogen 
(LN2) at various driving pressures with saturated and sub-cooled fluid. This example 
problem models one of the tests with an inlet boundary of saturated LH2 at 74.97 psia 
and -411 °F. 

Figure 14 shows the GFSSP model of the chilldown experiment. The 200 ft pipe line has 
been discretized into 30 pipe segments, each with 80 inches of pipe length. There are 31 
fluid nodes, each fluid node is connected with solid nodes. The total mass is distributed to 
3 1 solid nodes. The fluid and solid nodes are connected by solid to fluid conductors. The 
solid nodes are connected by solid to solid conductors. The boundary nodes 1 and 33 
represent the inlet dewar and ambient outlet respectively. 

The solid nodes are connected to the fluid nodes by Fluid-to-Solid Conductors, which 
model convection from the fluid to the pipe wall. The built-in Miropolskii correlation [17] 
is used to calculate the convection coefficient for the two-phase flow. Because the pipe is 
vacuum-jacketed, heat transfer between the pipe walls and the ambient is assumed 
negligible. 
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Figure 13. NBS Test Set-up of Cryogenic Transfer Line 
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Figure 14. GFSSP Model of the Cryogenic Pipeline 

Figure 15 shows the comparison between predicted and measured temperature at 4 
measuring stations are shown in this slide. The predictions and measurements compare 
reasonably well. In general measurements show more rapid chilldown than predictions. 
This discrepancy can be attributed to the deficiency of heat transfer correlation that does 
not account for other boiling regimes such as nucleate and transition. 

The comparison of measured and predicted chilldown time of the transfer line at different 
driving pressures for both liquid hydrogen and nitrogen is shown Table 6. Chilldown time 
decreases with increasing pressure primarily due to higher flowrate at higher pressure. 
Sub-cooling helps reduce chilldown time. Generally, predicted chilldown time is slightly 
higher than measured data. Discrepancy between the prediction and measurements can be 
attributed to the inaccuracy in heat transfer coefficient correlation as discussed in the 
previous paragraph. More details about this validation appear in Reference 18. 
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Figure 15. GFSSP’s Predicted Solid Temperatures (°F) Compared to Measurements 


Table 6 Measured and Predicted Chilldown Time for NBS Test Setup 


Saturated LH 2 chilldown time for various driving 
pressures 


Driving 

Pressure 

(paia) 

Saturation 

Temperature 

m 

Experimental 

Chilldown 

Time 

w 

Predicted 

Chilldown 

Ttme 

<•) 

74.97 

-411.06 

68 

70 

86.73 

-409.08 

62 

69 

111.72 

-406.4 

42 

50 

161.72 

-402.13 

30 

33 


Subcooled LH 2 chilldown time for various driving 
pressures. LH 2 is subcooled at -424.57 °F 


Driving 

Pressure 

ipsia) 

Experimental 

Chilldown 

Time 

{§) 

Predicted 

Chilldown 

Time 

<*) 

36 75 

148 

150 

61,74 

75 

80 

86.73 

62 

60 

111.72 

41 

45 

136.72 

32 

35 

161.7 

28 

30 


Saturated LN 2 chilldown time for various driving 
pressures 


Driving 

Pressure 

(psia) 

Saturation 

Temperature 

£Q 

Experimental 

Chilldown 

Time 

ja 

Predicted 

Chilldown 

Time 

(») 

8174 

-294,09 

165 

185 

74.97 

-289,71 

150 

160 

66.73 

-286.24 

130 

140 


Subcooled LN 2 chilldown time for various driving 
pressures. LN 2 is subcooled at -322.87 °F 


Driving 

Pressure 

ipsia) 

Experimental 

Chilldown 

Time 

(») 

Predicted 

Chilldown 

Time 

«•) 

36,75 

222 

250 

49,97 

170 

175 

61 74 

129 

140 

74.97 

100 

too 

86.73 

85 

90 
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4.4 Self-pressurization of cryogenic tank 

The purpose of this example is to demonstrate the simulation of self-pressurization of a 
Liquid Hydrogen Tank performed under the Multi-Purpose Hydrogen Test Bed (MHTB) 
program [19]. The purpose of the MHTB program is to test a Thermodynamic Vent 
System (TVS) to reduce boil-off in a Cryogenic Propellant Tank for long term storage of 
propellant in space as shown schematically in Figure 245. 



Figure 16. Thermodynamic Vent System in MHTB Ta nk 

The MHTB 5083 aluminum tank is cylindrical in shape with a height and diameter of 10 
feet and elliptic domes in both ends as shown in Figure 17. It has an internal volume of 639 
ft 3 and surface area of 379 ft 2 . Initially the tank is allowed to self- pressurize due to boil-off 
and by not allowing the vapor to vent. Once the pressure reaches the maximum allowable 
pressure, liquid hydrogen is introduced into the tank through the spray bar. The pressure 
starts falling due to heat transfer, and when the pressure reaches the minimum allowable 
pressure, the spray is stopped and the tank is allowed to self-pressurize and thus TVS cycle 
continues. The purpose of the GFSSP model is to simulate the initial self-pressurization 
when ullage pressure rises from the initial tank pressure to the upper bound pressure when 
the spray starts. The GFSSP model results were then compared with the test data. A 50% 
Fill Level case was modeled to simulate the self-pressurization test. 
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Figure 17. MHTB Test Tank and Supporting Hardware Schematic 





Figure 18. GFSSP Model of MHTB Test Tank 
Figure 18 shows the GFSSP model of self-pressurization in the MHTB Tank at the 50% fill 
level. Node 4 represents liquid hydrogen; Nodes 2, 8, 9, 10, and 1 1 represent the ullage at 
different fill levels. Node 3 is a pseudo-boundary node separating liquid hydrogen from 
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vapor hydrogen in the ullage space. Branches 45, 164, 162, 168, 169, 1610 and 1611 are for 
introducing liquid hydrogen into the tank through the TVS spray bar. These branches are 
currently inactive during self-pressurization of the tank. Node 7, 6, 12, 13, 14 and 15 are 
solid nodes representing the aluminum tank wall. Solid Node 7 is connected with liquid 
hydrogen stored in Fluid Node 4. In this model, heat leak through insulation is calculated 
in the User Subroutine and applied in the solid nodes as a source term. 

In this model a User Subroutine was used a) to model evaporative mass transfer at the 
liquid- vapor interface, b) to calculate the heat transfer coefficient between the wall and the 
fluid nodes, and c) to calculate heat transfer through the MLI blankets. The details of the 
modeling appear in reference 20. 

Figure 19 shows the comparison between GFSSP predictions (in green and blue) and the 
MHTB Test Data (in orange). GFSSP predictions of pressure are shown for a 
Degradation Factor of 1 and 2.8. The Degradation Factor is used to multiply Equation 
3.5.7 to represent the degradation of performance of the MLI. It is observed that a 
Degradation Factor of 2.8 matches the test data well. 
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Figure 19. Application Results for MHTB Self Pressurization Model 
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4.5 No vent chill and fill of cryogenic tank 

The purpose of this example is to demonstrate the simulation of the no-vent chill and fill 
method of chilling and filling a cryogenic tank. The practice of ta nk chilldown in micro- 
gravity environment is quite different than tank chilldown on the ground. On the ground, 
under normal gravity, a vent valve on top of the tank can be kept open to vent the vapor 
generated during the chilling process. The tank pressure can be kept close to atmospheric 
pressure while the tank is chilling down. In a micro-gravity environment, due to the 
absence of stratification, such practice may result in dumping a large amount of precious 
propellant overboard. The intent of the no- vent chill and fill method is to minimize the 
loss of propellant during chilldown of a propellant tank in a micro-gravity environment. 
The no-vent chill and fill method consists of a repeated cyclic process of charge, hold, 
and vent. 

During the charge cycle, a small quantity of liquid cryogen is injected into the evacuated 
tank. Some type of spray nozzle is usually used to break the incoming liquid into droplets. 
Initially, the liquid flashes due to the low tank pressure, and then the remaining liquid 
droplets evaporate as they contact warm hydrogen vapor or the tank wall. During the hold 
period, the circulating flow pattern induced from the spray nozzles provides convective 
heat transfer from cold vapor to the tank wall. The primary mode of heat transfer during 
the hold is convection. At the completion of the hold period, the pressure has risen 
considerably and the tank is ready to be vented. Since venting occurs as an isentropic 
blow-down, some additional cooling may be recovered by stage-wise venting. The key 
parameters of this method are (1) charge magnitude, (2) spray system selection, (3) mass 
flow rate, (4) hold duration, (5) acceleration environment, (6) desired tank wall 
temperature, and (7) maximum operating pressure. A reliable and inexpensive 
mathematical model will help designers to perform a large amount of calculations to 
optimize the key parameters. A GFSSP model was developed to simulate chilldown of the 
LH2 tank at the K-site Test Facility [21] and numerical predictions were compared with 
test data. 

The test set-up at the K-site Test Facility, shown in Figure 20, consists of a test tank, spray 
system, test tank valves, instrumentation, and the vacuum chamber. 
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Figure 20. K-Site Test Set-Up for No-Vent Fill Experiment 

The test tank selected was ellipsoidal with an 87 inch major diameter and a 1.2 to 1 major to 
minor axis ratio. The two ends are joined by a short 1.5 inch cylindrical section. The tank is 
made of 2219 aluminum chemically milled to a nominal thickness of 0.087 inches. Thicker 
sections exist where they were required for manufacturing (mainly weld lands). The tank has 
a 28.35 inch access flange on the top. The tank weighs 329.25 pounds, and the tank’s volume 
is 175 ft 3 . The tank was originally designed for a maximum operating pressure of 80 psia. 
Prior to the start of testing the tank was re-qualified by pneumatic test for a maximum 
operating pressure of 50 psia. The tank is covered with a blanket of 34 layers of multi-layer 
insulation (MLI) made with double aluminized Mylar and silk net spacers, and is supported 
by 12 fiberglass epoxy struts. The test environment ambient temperature was uniform and 
maintained at 530R ± 1R by an electrically heated shroud located outside the tank and in side 
the vacuum chamber. 

A nine node tank model (Figure 21) was developed to model this experiment. The tank was 
discretized into nine nodes and eight branches. Each fluid node was connected to a 
corresponding solid node. The total flowrate was equally distributed into 9 branches with 
fixed flowrate option. 
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Figure 21. GFSSP Nine Node Tank Model 
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Figure 22. Specified Inlet Flowrate (green) and Predicted Pressure (Orange) History 


Inlet flowrate into the tank and predicted pressure are shown in Figure 22. There are five 
short pulses of liquid hydrogen flow into the tank. After each pulse, there is a period of 
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no flow into the tank. During this period, the tank holds the propellant with vent valve 
closed for the tank to reach thermal equilibrium with cold vapor. Pressure rises almost 
instantaneously because of liquid turning into vapor with large increase in specific 
volume. The pressure continues to increase at a lower rate due to heat transfer from wall 
during this hold period. The hold period is followed by venting causing pressure to drop 
rapidly. After shutting the vent valve, the next pulse of inlet flow occurs. After five 
pulses, a constant inlet flow was maintained to fill the tank. Once the continuous filling 
starts, pressure initially drops due to some condensation of vapor in the tank. Once the 
tank is getting nearly filled, pressure rises due to the compression of vapor in the ullage 
space. 

Figure 22 shows the predicted mass history of hydrogen during the operation. There is 
very little hydrogen during the chilling process because of venting. The total amount of 
propellant vapor vented during this period is 32.5 lbs. This number compares well 
(within 1.5%) with measured propellant loss during the test. 

-H- [em 1 ]+ [em2] +[em^l+[em4]+ [em&]+ [em6]+[e m7]+[em8]+ [em9] 


700 



Figure 22. Predicted Hydrogen Mass History in the Tank 


Temperatures were measured during the test. The upper and lower bound of measured 
temperature history is shown in Figure 23. GFSSP model results for single and nine node 
are also shown to compare the predicted and test data. For nine node model, the 
centerline temperature is plotted. 
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GFSSP (Single Node) and GFSSP (9 Node Centerline) Wall Temperature 
Results Comparison to Test Wall Temperature Results 



Figure 23. Comparison of predicted and measured wall temperature during no- vent fill 

for the K-site test tank. 


5. Conclusions 

This paper explains the basic features of NASA’s Generalized Fluid System Simulation 
Program (GFSSP) and describes the additional capabilities of Version 6. Several 
numerical models are presented to illustrate code’s application in simulating Fluid 
Transient and Cryogenic Fluid Management (CFM) applications such as chilldown of 
cryogenic transfer line, no- vent chill and fill of cryogenic tank and self-pressurization of 
cryogenic tank. Numerical results are compared with test data. 
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